Quantitatively mapping local quality of super-resolution microscopy by rolling Fourier ring correlation

In fluorescence microscopy, computational algorithms have been developed to suppress noise, enhance contrast, and even enable super-resolution (SR). However, the local quality of the images may vary on multiple scales, and these differences can lead to misconceptions. Current mapping methods fail to finely estimate the local quality, challenging to associate the SR scale content. Here, we develop a rolling Fourier ring correlation (rFRC) method to evaluate the reconstruction uncertainties down to SR scale. To visually pinpoint regions with low reliability, a filtered rFRC is combined with a modified resolution-scaled error map (RSM), offering a comprehensive and concise map for further examination. We demonstrate their performances on various SR imaging modalities, and the resulting quantitative maps enable better SR images integrated from different reconstructions. Overall, we expect that our framework can become a routinely used tool for biologists in assessing their image datasets in general and inspire further advances in the rapidly developing field of computational imaging.

) of the background thresholding methods.During the rolling operation of the rFRC mapping, the intensity of center pixels from each block is summed (blue summation sign).The FRC value is calculated and assigned only if this summed value of the center pixels is larger than (blue tick sign) the threshold ('∑ > T'); otherwise, the center pixel is set to zero (blue cross sign) ('∑ < T').In this work, we provided two strategies for threshold determination.One is the user-defined hard threshold for the entire image ('15' as in this representative example).The other is the iterative wavelet transform method (yellow box), which automatically estimates the local threshold values.(b) rFRC maps using different background thresholds (c.f., Fig. 2a).(c) A representative SRRF data (left) (c.f., Fig. S7c) for illustration of two strategies of background thresholding (middle for hard threshold and right for adaptive threshold).the high-dimensional space (left panel of Fig. S13a).Thus, recovering the hidden ground-truth signals from the measured data is an ill-posed problem, while each reconstruction represents an approximation deviating from the real-world object to some extent (right panel of Fig. S13a).If two measurements are conducted on the same sample and considering the proper reconstruction model used, the distance between the two reconstructions may correlate to the data uncertainty (unbiased estimation, Fig. S13b).In this sense, we can use this distance (Reconstruction 1 versus Reconstruction 2) to approximate the magnitude of the reconstructed error.
To test this hypothesis, we restrict the derivation under specific and common circumstances for fluorescence imaging by a wide-field (WF) microscope, where the corresponding image can be processed by the minimalist form of Wiener deconvolution 9 .
Lemma.If we observe the object with an aberration-free, bandwidth-unlimited optical system, in which the object is sampled with an infinite sampling rate without noise, the minimalist form of Wiener deconvolution can recover the object perfectly: where ω and ω are the spatial and Fourier space coordinates, and I, o, h, and H represent the illumination, object, point spreading function (PSF), and its Fourier transform of the microscope.F andF  denote the Fourier and inverse Fourier transform operators, respectively.However, when considering the combined effects of sampling and noise, the model should be expressed as: where S and N denote the sampling and noise model, respectively, and imageWF is the final image collected by the WF microscope.This imaging model is considered closer to the status in the real physical world, as shown in Eq. ( 2), and it is far different from that of the minimalist Wiener model.

  ( ) ( ) ( )
As in Eq. ( 1), we still use the minimalist form of Wiener deconvolution with an aberration-free PSF (unbiased estimation) to process imageWF, which is given as: The corresponding result imagedecon in Eq. ( 4) is evidently at some distance from the Real object.This Lemma implies that the distance between the Real object and the reconstruction result is mainly caused by the combined effects of the sampling rates and mixture noise.
Corollary.Suppose we can control the variables to image the same object and capture a statistically independent image pair.In that case, the difference between the reconstruction and the Real object can be approximated by the distance between reconstructions from the image pair.In other words, this distance (Wiener result1 versus Wiener result2) can be related to the real errors of reconstruction (Wiener result1 and Wiener result2 versus Real object).
Then, the formula can be given as: where  denotes the union operation (defined as addition or multiplication in this work), and , Wiener Wiener D represents the distance between the Wiener result1 and Wiener result2.The Wienern is defined as: where 1 N , 2 N , 1 S and 2 S are two separate noise models ( 1 N and 2 N ) and sampling models ( 1 S and 2 S ).It is worth noting that the actual model of the union operation from Eq. ( 6) in real-world is complicated, and the combination of the two distances (reconstructions versus real objects) is challenging to express explicitly.
For simplicity, we used the Euclidean distance or multiplication in the cross-correlation distance as approximations.
Euclidean distance.In this example, we use the Euclidean distance to define the distance between two reconstructions and the addition operation to combine the distances.Then, for convenience, we modified Eq.
(5) by taking the Fourier transform of the Wiener results and Real object: ,   represents the Euclidean distance.The image is captured with an infinite sampling rate, and it is corrupted with additive noise only, denoted as nn: By calculating the Euclidean distance, we can simplify Eq. ( 7) to obtain the following form: where n1 and n2 denote the separate additive noise in the two observations.Then, we simplify Eq. ( 9) as: In Eq. ( 10) it clearly shows the close relationship between the distance between the two reconstructions (Wiener result1 versus Wiener result2) and the real errors of the two reconstructions (Wiener result1 and Wiener result2 versus Real object).
Cross-correlation distance.Alternatively, if we use the cross-correlation to calculate the distance, the addition operation in the Euclidean distance calculation needs to be changed to a multiplication operation: where F , and I F represent the Fourier transforms of Wiener result1, Wiener result2, and the Real object, respectively.Then, the following formula can be derived: where are defined.Similar to the Euclidean distance, Eq. ( 12) shows the close relationship between the cross-correlation (Wiener result1 versus Wiener result2) and the multiplication of the cross-correlation (Wiener result1 versus Real object and Wiener result2 versus Real object).
As shown above, both the Euclidean distance and cross-correlation can be used to estimate the distance between real objects and the corresponding reconstructions.However, it is not easy to choose a superior one to quantitatively map these distances in multiple dimensions (2D image or 3D volume).Conventional spatial methods, such as the spatial subtraction in the RSM 10 , calculate the 'absolute differences' 11 .These methods are prone to false negatives in the distance map upon intensity fluctuations and sample movements 11 .
Therefore, we introduced a method to measure the distance between two signals in the Fourier domain, namely, Fourier ring correlation (FRC, Methods) [12][13][14] or spectral signal-to-noise ratio (SSNR), describing the highest acceptable frequency component between two signals.The FRC metric has been used as a practical resolution criterion for super-resolution (SR) fluorescence and electron microscopy.Here, we used it to quantify the similarity (distance) between two signals for its insensitivity to intensity changes and micromovements, and also for its quantifying the 'relative error' or 'saliency-based error' (the highest credible frequency component).Because the aberrations in the system may not change during these two independent observations, these aberration-induced errors (biased estimation, Fig. S13b) are difficult to estimate when using simple spatial or frequency methods.Instead, these errors may be visible by the FRC since it defines the most reliable frequency component.Overall, the FRC can be a superior choice to quantify the distance between two signals.Furthermore, in this work, considering the FRC as a global similarity estimation between two images, we extended the FRC metric to a rolling FRC (rFRC, Methods) map to provide local distance measurements at the pixel level.Therefore, we can quantitatively map the errors in the multidimensional reconstructed signals without the ground truth.
We compared our developed Fourier domain method, rFRC mapping, with a spatial domain method, standard deviation (STD), for data uncertainty measurement in the application of Richardson-Lucy (RL) deconvolution (Fig. S14).We repetitively added independent noise 8 times and then applied RL deconvolution on the resulting 8 raw wide-field images.We first calculated the STD on the 2 deconvolved images.Because of severe intensity fluctuations, it is found that the STD with 2 images failed to highlight the three different noise levels.As a result, for this spatial domain method, at least 8 images for STD calculation are required to assess the data uncertainty without amplifying false negatives.By contrast, the rFRC maps using 2 and 8 images exhibit almost identical distributions and both successfully capture the three levels of reconstruction quality.This test suggests that our rFRC mapping method is more stable and efficient compared to the commonly used spatial domain method.
If the two signals lose the identical component in reconstructions, the rFRC may indicate a false positive.This is a potential pitfall of the rFRC method, which could be compensated for by the modified resolution scaled error map (RSM, Methods) method 10 .Specifically, the RSM is built upon three major assumptions: Because the RSM requires a diffraction-limited WF image as the reference, it must have a sufficiently high SNR.Otherwise, the noise contained in the WF image may induce false-negative effects.
(ii).The transformation from the SR scale to the low-resolution (LR) scale is approximated as a global spatially invariant Gaussian kernel convolution.If the Gaussian-shaped PSF is not spatially invariant, estimation errors will occur.
(iii).The illumination, emission fluorescence intensity, and background are homogeneous.The produced RSM is weighted by the background, illumination intensity, and emission fluorescence intensity of the corresponding label.
If the assumptions mentioned above do not follow, the RSM may introduce false negatives in the estimated error maps.Assumptions (i) and (ii) may be satisfied under the configurations of single-molecule localization microscopy (SMLM).For example, the WF reference image of SMLM is created by averaging thousands of blinking frames, which removes the noise and generates a high SNR reference.In addition, because the routine imaging field has a small field of view (FOV), it can be regarded as uniform illumination.
The RSM can only find errors on an LR scale (usually error components of large magnitudes), such as misrepresentation or disappearance of structures.In other words, errors in the SR scale (small-magnitude error components) highlighted by the RSM may not be correct.Therefore, we segment the normalized RSM (value less than 0.5 is set to 0) to remove the corresponding components of small magnitude (Methods), leaving possible errors of large magnitude (such as the complete absence of structure) to complement the rFRC map.
Finally, we merged the developed rFRC map and the segmented normalized RSM in the green and red channels to generate a complete PANEL visualization (Methods).the distributions of FRC are just separatable.Paired lines became more separable as overlaps decreased in 8pixel and 16-pixel cases and were distinct in 32-pixel and 64-pixel cases.Images must satisfy the Nyquist sampling criteria to achieve maximal resolution, so their point spread function (PSF) should cover at least 3pixel.Therefore, the separation of rFRC of paired lines 4-pixel apart means the minimum detectable scale of rFRC map is up to its limit.By involving the rolling operation, we have addressed a major limitation of the previous FRC map 10 , which is challenging to correlate the block-wise map to the SR image content (Fig. S15c, S15d, and S17).The FRC curve (black), 3σ threshold curve (blue), and 1/7 threshold curve (red) for a 64 × 64 pixels image (solid) and a 512 × 512 pixels image (dashed).For an image with a large size (512 × 512 pixels), the 1/7 threshold attains a similar result to the 3σ curve criterion (green circle).However, for a small image (64 × 64 pixels), the 1/7 threshold is smaller than all correlation values in the FRC curve, failing to yield the cutoff frequency.(b) The uncertainty of FRC calculation using different block sizes by 1/7 threshold curve (left) and 3σ threshold curve (right).We downsampled the 2D-STORM captured microtubule image (c.f., Fig. 3f, 10 nm pixel size, 4096 pixel-number) with 16, 32, 64, and 128 times to create different image sizes and convoluted the resulting images with a 120 nm PSF.After that, Poisson and 5% Gaussian noise were injected into the image.This procedure was repeated 20 times independently and the FRC calculations were performed with different criteria.(c) rFRC maps using different block sizes (c.f., Fig. 2a).Although the smaller block size (e.g., 32 × 32 pixels) may enable finer mapping, the overall distributions of these rFRC resolution maps using different block sizes are close to each other.On the other hand, the overly small block size may lead to an overconfident resolution value and larger uncertainty.Therefore, to balance the compromise between mapping scale and estimation stability, we chose a block size of 64 × 64 pixels as default in this work.(d) FRC maps using different block sizes.Scale bar: 500 nm.

Fig. S1 |
Fig. S1 | Abstract workflow.(a) Abstract workflow.Only when the corresponding tasks satisfy two conditions, i.e., (i) belonging to 2D data and (ii) the existence of a wide-field reference, will the RSM be included in the PANEL visualization.(b) Our framework for estimating different types of uncertainties.At the SR scale, our method is capable of mapping (i) data uncertainty of image reconstructions without referencing the ground-truth (Reconstruction-1 vs. Reconstruction-2); (ii) large scale errors induced by model-bias referencing wide-field image (Reconstruction vs. Wide-field reference); and (iii) full error of reconstructions/predictions with ground-truth (Reconstruction vs. Ground-truth).

Fig. S2 |
Fig. S2 | Two background skip strategies for rFRC mapping.(a) Workflows (c.f., Fig.2a) of the background thresholding methods.During the rolling operation of the rFRC mapping, the intensity of center pixels from each block is summed (blue summation sign).The FRC value is calculated and assigned only if this summed value of the center pixels is larger than (blue tick sign) the threshold ('∑ > T'); otherwise, the center pixel is set to zero (blue cross sign) ('∑ < T').In this work, we provided two strategies for threshold determination.One is the user-defined hard threshold for the entire image ('15' as in this representative example).The other is the iterative wavelet transform method (yellow box), which automatically estimates the local threshold values.(b) rFRC maps using different background thresholds (c.f., Fig.2a).(c) A representative SRRF data (left) (c.f., Fig.S7c) for illustration of two strategies of background thresholding (middle for hard threshold and right for adaptive threshold).

Fig. S3 |
Fig. S3 | Color maps for map display and Otsu threshold for PANEL pinpointing.(a) The representative color-coded images and color indexes of jet (left), black jet (middle), and shifted jet (right) color maps.The image is adapted from Fig. S7a.(b) Otsu threshold for PANEL highlighting.Left: The rFRC map of the SRRF dataset in Fig. S7c, displayed in the sJet color map.Middle and right: The Full rFRC map (middle) and the rFRC map after the Otsu threshold (right), regions with low reliability in the SRRF reconstruction are pointed by green.

Fig. S7 |
Fig. S7 | Open-source 2D-SMLM and SRRF experimental datasets evaluations.(a) From left to right: MLE localization result of 500 high-density images of tubulins from the EPFL website (Methods); the rFRC map of the MLE; full merged RSM and rFRC map of the MLE; PANEL visualization.(b) From left to right: Corresponding wide-field image; MLE image convolved back to its original low-resolution scale; RSM of the MLE; FRC map of the MLE.(c) From left to right: Diffraction-limited TIRF image; SRRF reconstruction result of 100 fluctuation images (GFP-tagged microtubules in live HeLa cells, Methods); rFRC map of SRRF; PANEL visualization.Scale bar: 2 μm.

Fig. S9 |
Fig. S9 | Another representative example of STORM fusion (COS-7 cells, heavy chain clathrin-coated pits labeled with Alexa Fluor 647).(a) The rFRC map of ME-MLE (top), the superiority map (middle) for fusion, and the TIRF image (bottom) are shown on the left.The rFRC maps of the SE-Gaussian (top) and fusion (middle) results, and the fusion result ('Fused', bottom) are displayed on the right.We found that the ME-MLE method achieves superior performance in the regions containing a strong background and that the SE-Gaussian method obtains better reconstruction quality in the regions containing a weak background.(b) Magnified results for a single CCP of ME-MLE (top), SE-Gaussian (middle), and fusion ('Fused,' bottom) are shown on the left, and the corresponding rFRC maps are demonstrated on the right.The mean resolutions are marked on the top left of the rFRC maps.In addition to the stable performance of fusion in the whole field of view, as highlighted in (a), the rFRC map assists in fusing fine structures such as a single ring-shaped CCP, enabling higher mean resolution.Scale bars: (a) 5 μm; (b) 100 nm.

Fig. S13 |
Fig. S13 | rFRC distance.(a) The basic concept of rFRC evaluation.The observer in the original domain indicates the captured raw images, and the reconstruction in the target domain represents the reconstructed SR images.(b) The data error, the model error, and the rFRC distance.The unbiased estimation indicates that the reconstruction model accurately describes the corresponding real-world model.Then the expected value of reconstructions will be identical to the actual value (green circle).The biased estimation represents that the reconstruction model is different from the real-world one.The expected value (purple circle) of reconstructions will deviate from the actual value (black arrow, model error).

Fig. S14 |
Fig. S14 | Evaluations in spatial domain versus Fourier domain.A series of filaments were convoluted with a wide-field PSF (numerical aperture as 1.4).We gradually decreased the noise level from left to right in three levels (left bottom).Left top: RL deconvolution result.Middle: STD (standard deviation) results of 2 (top) and 8 (bottom) RL deconvolved images.Middle: rFRC map of 2 (top) and 8 (bottom) RL deconvolved images.The 8 images were splitted as 4 batches for rFRC mapping, and the 4 resulting rFRC maps were averaged to create the final rFRC map (8).Scale bar: 500 nm.

Fig. S15 |
Fig. S15 | The stability of rFRC map.(a)The FRC curve (black), 3σ threshold curve (blue), and 1/7 threshold curve (red) for a 64 × 64 pixels image (solid) and a 512 × 512 pixels image (dashed).For an image with a large size (512 × 512 pixels), the 1/7 threshold attains a similar result to the 3σ curve criterion (green circle).However, for a small image (64 × 64 pixels), the 1/7 threshold is smaller than all correlation values in the FRC curve, failing to yield the cutoff frequency.(b) The uncertainty of FRC calculation using different block sizes by 1/7 threshold curve (left) and 3σ threshold curve (right).We downsampled the 2D-STORM captured microtubule image (c.f., Fig.3f, 10 nm pixel size, 4096 pixel-number) with 16, 32, 64, and 128 times to create different image sizes and convoluted the resulting images with a 120 nm PSF.After that, Poisson and 5% Gaussian noise were injected into the image.This procedure was repeated 20 times independently and the FRC calculations were performed with different criteria.(c) rFRC maps using different block sizes (c.f., Fig.2a).Although the smaller block size (e.g., 32 × 32 pixels) may enable finer mapping, the overall distributions of these rFRC resolution maps using different block sizes are close to each other.On the other hand, the overly small block size may lead to an overconfident resolution value and larger uncertainty.Therefore, to balance the compromise between mapping scale and estimation stability, we chose a block size of 64 × 64 pixels as default in this work.(d) FRC maps using different block sizes.Scale bar: 500 nm.

Fig. S16 |
Fig.S16| The resolvability of rFRC map.We simulated structures that contained pairs of lines with spacing gradually increases, i.e., 2, 4, 8, 16, 32, 64 pixels (pointed by red arrows), and convoluted them by a PSF with a 4-pixel FWHM (pixel size 10 nm).To test the maximum resolvability of rFRC, we included different noise levels on the two lines.Specifically, we added 10% and 50% Gaussian noise on the left and right lines, respectively.After that, we applied rFRC mapping on the resulting images (left panel) and calculated the FRC value distributions (right panel) of pixels on the left (yellow) and right lines (green).